library(cjoint);library(dplyr)
library(FindIt)
setwd("~/Dropbox/Gender Quota in Korea/_Pink-CollarWorker/Data/UK_analysis")
library(foreign)
library(readstata13)
library(tidyverse)
library(cowplot)

#### Figure 3: All Respondents AMCE
## (1) Candidate Support ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.3", "Q6.3", "Q7.3"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = factor(Education, levels = c("High school graduate", "University graduate", "Post-graduate (Master or Doctoral degree)")),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))

model <- amce(selected ~  Age + Education + Working + Gender + Years.in.Politics, data = data, cluster=TRUE, respondent.id = "Response.ID")
model$estimates

labels <- c("Gender", "    Male", "    Female",
            "Occupational Background", "    White Collar","    Working Class",
            "Age", "   30", "   40", "   50", "   60", 
            "Education", "    High school graduate", "    University graduate", "    Post-graduate (Master or Doctoral degree)", 
            "Years in Politics", "    0", "    1", "    3", "    8"
)


estimates <- c(NA, 0, model$estimates[3][[1]][[1]], # gender
               NA, 0, model$estimates[4][[1]][[1]], # class
               NA, 0, model$estimates[1][[1]][[1]], model$estimates[1][[1]][[3]],model$estimates[1][[1]][[5]],
               NA, 0, model$estimates[2][[1]][[1]], model$estimates[2][[1]][[3]],
               NA, 0, model$estimates[5][[1]][[1]],  model$estimates[5][[1]][[3]],  model$estimates[5][[1]][[5]]
)

se <- c(NA, 0, model$estimates[3][[1]][[2]],
        NA, 0, model$estimates[4][[1]][[2]],
        NA, 0, model$estimates[1][[1]][[2]], model$estimates[1][[1]][[4]],model$estimates[1][[1]][[6]],
        NA, 0, model$estimates[2][[1]][[2]], model$estimates[2][[1]][[4]],
        NA, 0, model$estimates[5][[1]][[2]],  model$estimates[5][[1]][[4]],  model$estimates[5][[1]][[6]]
)

face <- c("bold", "plain", "plain",
          "bold", "plain", "plain",
          "bold", "plain", "plain", "plain", "plain",
          "bold", "plain", "plain", "plain",
          "bold", "plain", "plain", "plain", "plain")

type = c("gender", "gender", "gender",
         "occupation", "occupation", "occupation", 
         "age", "age","age","age", "age",
         "education", "education", "education", "education",
         "year", "year", "year", "year", "year")

support <- data.frame(labels = factor(labels, levels = rev(labels)), 
                      estimates = estimates,
                      se = se,
                      face = face,
                      type = type)

support <- support %>% mutate(ul = estimates + 1.96 * se,
                              ll = estimates - 1.96 * se,
                              ref = ifelse(estimates == 0 | is.na(estimates), "reference", "normal"),
                              dv = "Candidate Preference")

support
## (2) Competency: Advocate R's Priority Issue ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.7", "Q6.7", "Q7.7"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Female"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))
model <- amce(selected ~  Age + Education + Working + Gender + Years.in.Politics, data = data, cluster=TRUE, respondent.id = "Response.ID")


labels <- c("Gender", "    Male", "    Female",
            "Occupational Background", "    White Collar","    Working Class",
            "Age", "   30", "   40", "   50", "   60", 
            "Education", "    High school graduate", "    University graduate", "    Post-graduate (Master or Doctoral degree)", 
            "Years in Politics", "    0", "    1", "    3", "    8"
)


estimates <- c(NA, 0, model$estimates[3][[1]][[1]], # gender
               NA, 0, model$estimates[4][[1]][[1]], # class
               NA, 0, model$estimates[1][[1]][[1]], model$estimates[1][[1]][[3]],model$estimates[1][[1]][[5]],
               NA, 0, model$estimates[2][[1]][[1]], model$estimates[2][[1]][[3]],
               NA, 0, model$estimates[5][[1]][[1]],  model$estimates[5][[1]][[3]],  model$estimates[5][[1]][[5]]
)

se <- c(NA, 0, model$estimates[3][[1]][[2]],
        NA, 0, model$estimates[4][[1]][[2]],
        NA, 0, model$estimates[1][[1]][[2]], model$estimates[1][[1]][[4]],model$estimates[1][[1]][[6]],
        NA, 0, model$estimates[2][[1]][[2]], model$estimates[2][[1]][[4]],
        NA, 0, model$estimates[5][[1]][[2]],  model$estimates[5][[1]][[4]],  model$estimates[5][[1]][[6]]
)

face <- c("bold", "plain", "plain",
          "bold", "plain", "plain",
          "bold", "plain", "plain", "plain", "plain",
          "bold", "plain", "plain", "plain",
          "bold", "plain", "plain", "plain", "plain")

type = c("gender", "gender", "gender",
         "occupation", "occupation", "occupation", 
         "age", "age","age","age", "age",
         "education", "education", "education", "education",
         "year", "year", "year", "year", "year")

issue <- data.frame(labels = factor(labels, levels = rev(labels)), 
                    estimates = estimates,
                    se = se,
                    face = face,
                    type = type)

issue <- issue %>% mutate(ul = estimates + 1.96 * se,
                          ll = estimates - 1.96 * se,
                          ref = ifelse(estimates == 0 | is.na(estimates), "reference", "normal"),
                          dv = "Competency")

results <- issue %>%
  bind_rows(support)

results <- results %>%
  mutate(dv = fct_relevel(dv, "Competency", "Candidate Preference"))
results

ggplot(results, aes(x=labels, y=estimates, ymin=ll, ymax=ul, shape= type)) +
  labs(y = "Effect on Pr(Being Selected)") +
  geom_pointrange() + 
  theme_bw()  +
  coord_flip() + 
  scale_x_discrete(name="",labels=rev(labels)) +
  facet_wrap(vars(dv), nrow = 1) +
  geom_hline(yintercept = 0,colour="darkgrey",linetype=1) +
  #theme(axis.text.y = element_text(hjust=0)) +
  theme(
    axis.text.y = element_text(colour = "black",
                               hjust = 0 , vjust=.5,
                               face = rev(face), size = 14),
    axis.ticks.y = element_blank(),
    strip.text = element_text(size = 14),
    legend.position = "none"
  ) + scale_shape_manual(values = c(18, 2, 20, 15, 7 ))

## Figure 4: Average Marginal Interaction Effect 
#### Candidate Support (All) ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.3", "Q6.3", "Q7.3"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = factor(Education, levels = c("High school graduate", "University graduate", "Post-graduate (Master or Doctoral degree)")),
                        # Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))

data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit1 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID,
                    pair.id = data$pairid, diff = TRUE, nway=2)
ConditionalEffect(fit1, treat.fac="Working", cond.fac="Gender")

labels <- c( "All", "All")

cond.estimates<- ConditionalEffect(fit1, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
cond.estimates
estimates <- c(cond.estimates[[1]][2,1],cond.estimates[[2]][2,1])
estimates
ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])
condition <- c("Candidate Gender = Male",  "Candidate Gender = Female")
type = c("male",  "female")

support.all <- data.frame(labels = factor(labels), 
                      estimates = estimates,
                      ll=ll,
                      ul=ul,
                      type=type, 
                      condition=condition)

support.all <- support.all %>% mutate(
                              dv = "Candidate Support")

support.all

#### Candidate Support (Labour) ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.3", "Q6.3", "Q7.3"),
                       covariates = c("Q3.1"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% filter(Q3.1=="2")
data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = factor(Education, levels = c("High school graduate", "University graduate", "Post-graduate (Master or Doctoral degree)")),
                        # Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))

data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit2 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID, nway=2,
                    pair.id = data$pairid, diff = TRUE)
ConditionalEffect(fit2, treat.fac="Working", cond.fac="Gender")

labels <- c( "Labour",
             "Labour")

cond.estimates<- ConditionalEffect(fit1, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
cond.estimates
estimates <- c(cond.estimates[[1]][2,1], cond.estimates[[2]][2,1])

ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])

condition <- c("Candidate Gender = Male", 
               "Candidate Gender = Female")
type = c("male", 
         "female")

support.lab <- data.frame(labels = factor(labels), 
                          estimates = estimates,
                          ll=ll,
                          ul=ul,
                          type=type, 
                          condition=condition)
support.lab <- support.lab %>% mutate(dv = "Candidate Preference")

## Conservative ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.3", "Q6.3", "Q7.3"),
                       covariates = c("Q3.1"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% filter(Q3.1=="1")
data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = factor(Education, levels = c("High school graduate", "University graduate", "Post-graduate (Master or Doctoral degree)")),
                        # Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))

data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit3 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID, nway=2,
                    pair.id = data$pairid, diff = TRUE)
ConditionalEffect(fit3, treat.fac="Working", cond.fac="Gender")

labels <- c( "Conservative", "Conservative")

cond.estimates<- ConditionalEffect(fit3, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
estimates <- c(cond.estimates[[1]][2,1], cond.estimates[[2]][2,1])

ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])

condition <- c("Candidate Gender = Male", "Candidate Gender = Female")
type = c("male", "female")

support.con <- data.frame(labels = factor(labels), 
                          estimates = estimates,
                          ll=ll,
                          ul=ul,
                          type=type, 
                          condition=condition)

support.con <- support.con %>% mutate(dv = "Candidate Preference")

support.con$party <- "Conservative"
support.lab$party <- "Labour"
support.all$party <- "All"
support <- support.all %>% bind_rows(support.lab) %>% bind_rows(support.con)

support <- support %>% mutate(labels = fct_relevel(labels, "Conservative", "Labour", "All"))
p1 <- ggplot(support, aes(x=labels, y=estimates, ymin=ll, ymax=ul, shape = party)) +
  labs(y = "Effect on Pr(Being Selected)") +
  geom_pointrange() + 
  theme_bw()  +
  coord_flip(ylim = c(-0.2, 0.2)) + 
  scale_x_discrete(name="") +
  facet_wrap(vars(condition), nrow = 2) +
  geom_hline(yintercept = 0,colour="darkgrey",linetype=1) +
  theme(
    axis.text.y = element_text(colour = "black",
                               hjust = 0 , vjust=.5,
                               size = 14),
    axis.ticks.y = element_blank(),
    strip.text = element_text(size = 14),
    legend.position = "none"
  )   + ggtitle("Working-Class Effect on Preference")
p1


#### Competency  (All) ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.7", "Q6.7", "Q7.7"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = factor(Education, levels = c("High school graduate", "University graduate", "Post-graduate (Master or Doctoral degree)")),
                        # Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))

data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit1 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID,
                    pair.id = data$pairid, diff = TRUE, nway=2)
ConditionalEffect(fit1, treat.fac="Working", cond.fac="Gender")

labels <- c( "All", "All")

cond.estimates<- ConditionalEffect(fit1, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
cond.estimates
estimates <- c(cond.estimates[[1]][2,1],cond.estimates[[2]][2,1])
estimates
ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])
condition <- c("Candidate Gender = Male",  "Candidate Gender = Female")
type = c("male",  "female")

issue.all <- data.frame(labels = factor(labels), 
                          estimates = estimates,
                          ll=ll,
                          ul=ul,
                          type=type, 
                          condition=condition)

issue.all <- issue.all %>% mutate(
  dv = "Competency")


## Competency- Labour ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.7", "Q6.7", "Q7.7"),
                       covariates = c("Q3.1"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% filter(Q3.1=="2")
data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))
data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit2 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID, nway=2,
                    pair.id = data$pairid, diff = TRUE)

cond.estimates<- ConditionalEffect(fit2, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
labels <- c( "Labour", "Labour")

estimates <- c(cond.estimates[[1]][2,1] 
               ,cond.estimates[[2]][2,1])

ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])
condition <- c("Candidate Gender = Male", "Candidate Gender = Female")
type = c("male","female")

issue.lab <- data.frame(labels = factor(labels), 
                        estimates = estimates,
                        ll=ll,
                        ul=ul,
                        type=type, 
                        condition=condition)

issue.lab <- issue.lab %>% mutate(
  dv = "Competency")

issue.lab$party <- "Labour"

#### Competency - Conservative ####
data <- read.qualtrics("Working-Class+Women+Conjoint+-+UK_February+9,+2022_03.16.csv",
                       responses = c("Q5.7", "Q6.7", "Q7.7"),
                       covariates = c("Q3.1"),
                       respondentID = "ResponseId", new.format = TRUE)

data <- data %>% filter(Q3.1=="1")
data <- data %>% drop_na(selected)
data <- data %>% mutate(Education = fct_relevel(Education, "High school graduate"),
                        Gender = fct_relevel(Gender, "Male"),
                        Working = factor(case_when(Former.Occupation == "Janitor" ~ "Working Class",
                                                   Former.Occupation == "Retail clerk" ~ "Working Class",
                                                   Former.Occupation == "Lawyer" ~ "White Collar",
                                                   Former.Occupation == "Business entrepreneur" ~ "White Collar")),
                        Working = fct_relevel(Working, "White Collar"))
data <- data %>% mutate(pairid = paste(respondent, task, sep="_"))

fit3 <- CausalANOVA(formula = selected ~  Age + Education + Working + Years.in.Politics + Gender , 
                    int2.formula = ~Working:Gender, data=data,  cluster = data$Response.ID, nway=2,
                    pair.id = data$pairid, diff = TRUE)

cond.estimates<- ConditionalEffect(fit3, treat.fac="Working", cond.fac="Gender")$ConditionalEffects
labels <- c( "Conservative", "Conservative")



estimates <- c(cond.estimates[[1]][2,1] 
               ,cond.estimates[[2]][2,1])

ll <- c(cond.estimates[[1]][2,3],cond.estimates[[2]][2,3])
ul <- c(cond.estimates[[1]][2,4],cond.estimates[[2]][2,4])
condition <- c("Candidate Gender = Male", "Candidate Gender = Female")
type = c("male","female")

issue.con <- data.frame(labels = factor(labels), 
                        estimates = estimates,
                        ll=ll,
                        ul=ul,
                        type=type, 
                        condition=condition)

issue.con <- issue.con %>% mutate(
  dv = "Competency")
issue.con$party <- "Conservative"
issue.all$party <- "All"
issue <- issue.all %>% bind_rows(issue.lab) %>% bind_rows(issue.con)

issue <- issue %>% mutate(labels = fct_relevel(labels, "Conservative", "Labour", "All"))

p2 <- ggplot(issue, aes(x=labels, y=estimates, ymin=ll, ymax=ul, shape=party)) +
 # scale_shape_manual(values = c(19,1)) +
  labs(y = "Effect on Pr(Being Selected)") +
  geom_pointrange() + 
  theme_bw()  +
  coord_flip(ylim = c(-0.2, 0.2)) + 
  scale_x_discrete(name="") +
  facet_wrap(vars(condition), nrow = 2) +
  geom_hline(yintercept = 0,colour="darkgrey",linetype=1) +
  theme(
    axis.text.y = element_text(colour = "black",
                               hjust = 0 , vjust=.5,
                               size = 14),
    axis.ticks.y = element_blank(),
    strip.text = element_text(size = 14),
    legend.position = "none"
  ) +  ggtitle("Working-Class Effect on Competency")
p2
library(ggpubr)
ggarrange(p2, p1, ncol=2)


